library(rhdf5)
library(raster)
## Loading required package: sp
library(rgdal)
## rgdal: version: 1.1-10, (SVN revision 622)
## Geospatial Data Abstraction Library extensions to R successfully loaded
## Loaded GDAL runtime: GDAL 2.0.1, released 2015/09/15
## Path to GDAL shared files: C:/Users/sjgraves/Documents/R/win-library/3.3/rgdal/gdal
## Loaded PROJ.4 runtime: Rel. 4.9.2, 08 September 2015, [PJ_VERSION: 492]
## Path to PROJ.4 shared files: C:/Users/sjgraves/Documents/R/win-library/3.3/rgdal/proj
## Linking to sp version: 1.2-3
# set file path
# use ../ to move up a directory to the data folder
h5file <- "../NEONdata/D17-California/TEAK/2013/spectrometer/reflectance/Subset3NIS1_20130614_100459_atmcor.h5"
# view h5 structure
h5ls(h5file)
## group name otype dclass
## 0 / ATCOR_Input_File H5I_DATASET STRING
## 1 / ATCOR_Processing_Log H5I_DATASET STRING
## 2 / Aerosol Optical Depth H5I_DATASET INTEGER
## 3 / Aspect H5I_DATASET FLOAT
## 4 / Cast Shadow H5I_DATASET INTEGER
## 5 / Dark Dense Vegetation Classification H5I_DATASET INTEGER
## 6 / Haze-Cloud-Water Map H5I_DATASET INTEGER
## 7 / Illumination Factor H5I_DATASET INTEGER
## 8 / Path Length H5I_DATASET FLOAT
## 9 / Processing Version H5I_DATASET STRING
## 10 / Reflectance H5I_DATASET INTEGER
## 11 / Shadow_Processing_Log H5I_DATASET STRING
## 12 / Sky View Factor H5I_DATASET INTEGER
## 13 / Skyview_Processing_Log H5I_DATASET STRING
## 14 / Slope H5I_DATASET FLOAT
## 15 / Slope_Aspect_Processing_Log H5I_DATASET STRING
## 16 / Solar Azimuth Angle H5I_DATASET FLOAT
## 17 / Solar Zenith Angle H5I_DATASET FLOAT
## 18 / Surface Elevation H5I_DATASET FLOAT
## 19 / Visibility Index Map H5I_DATASET INTEGER
## 20 / Water Vapor Column H5I_DATASET INTEGER
## 21 / coordinate system string H5I_DATASET STRING
## 22 / flightAltitude H5I_DATASET FLOAT
## 23 / flightHeading H5I_DATASET FLOAT
## 24 / flightTime H5I_DATASET FLOAT
## 25 / fwhm H5I_DATASET FLOAT
## 26 / map info H5I_DATASET STRING
## 27 / to-sensor azimuth angle H5I_DATASET FLOAT
## 28 / to-sensor zenith angle H5I_DATASET FLOAT
## 29 / wavelength H5I_DATASET FLOAT
## dim
## 0 1
## 1 1
## 2 544 x 578 x 1
## 3 544 x 578 x 1
## 4 544 x 578 x 1
## 5 544 x 578 x 1
## 6 544 x 578 x 1
## 7 544 x 578 x 1
## 8 544 x 578 x 1
## 9 1
## 10 544 x 578 x 426
## 11 1
## 12 544 x 578 x 1
## 13 1
## 14 544 x 578 x 1
## 15 1
## 16 1 x 1
## 17 1 x 1
## 18 544 x 578 x 1
## 19 544 x 578 x 1
## 20 544 x 578 x 1
## 21 1
## 22 5732053
## 23 5732053
## 24 5732053
## 25 426 x 1
## 26 1
## 27 544 x 578 x 1
## 28 544 x 578 x 1
## 29 426 x 1
# import spatial information
mapInfo <- h5read(h5file,
"map info",
read.attributes = T)
# set reflectance info
reflInfo <- h5readAttributes(h5file,
"Reflectance")
# set scale factor and no data value objects
scaleFactor <- reflInfo$`Scale Factor`
noDataValue <- as.numeric(reflInfo$`data ignore value`)
NEON wrote functions to do many data extractions manually. Information is found on the website
rm(list = ls())
source("neon_hdf5_functions.R")
f <- "../NEONdata/D17-California/TEAK/2013/spectrometer/reflectance/Subset3NIS1_20130614_100459_atmcor.h5"
# define the CRS definition by EPSG code
epsg <- 32611
# set wavelengths
wavelengths <- h5read(f,"wavelength")
### final Code ####
# H5close()
# find the dimensions of the data to help determine the slice range
# returns cols, rows, wavelengths
dims <- get_data_dims(fileName = f)
# open band, return cleaned and scaled raster
band <- open_band(fileName=f,
bandNum = 56,
epsg=epsg)
# plot data
plot(band,
main="Raster for Lower Teakettle - B56")
# extract 3 bands
# create alist of the bands
bands <- list(58, 34, 19)
# use lapply to run the band function across all three of the bands
rgb_rast <- lapply(bands, open_band,
fileName=f,
epsg=epsg)
# create a raster stack from the output
rgb_rast <- stack(rgb_rast)
# plot the output, use a linear stretch to make it look nice
plotRGB(rgb_rast,
stretch='lin')
# CIR create alist of the bands
bands <- c(90, 34, 19)
CIRStack <- create_stack(f,
bands,
epsg)
plot_stack(CIRStack,
title="Color Infrared (CIR) Image")
# create a list of the bands
bands <- list(152,90,58)
aStack <- create_stack(f, bands, epsg)
plot_stack(aStack,
title="another combo")
# export as a GeoTIFF
writeRaster(CIRStack,
file="../outputs/TEAK/TEAK_CIR_image.tif",
format="GTiff",
overwrite=TRUE)